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D I Abstract. The dynamical behaviours of a kinetically constrained spin model 

H ■ (Fredrickson-Andersen model) on a Bethe lattice are investigated by a perturbation 

analysis that provides exact final states above the nonergodic transition point. It is 
^ ' observed that the time-dependent solutions of the derived dynamical systems obtained 

"^ ■ by the perturbation analysis become systematically closer to the results obtained by 

Monte Carlo simulations as the order of a perturbation series is increased. This 
systematic perturbation analysis also clarifies the existence of a dynamical scaling 
law, which provides a implication for a universal relation between a size scale and a 
' y ' time scale near the nonergodic transition. 
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Systematic perturbation approach for a dynamical scaling law in a kinetically constrained spin modeU 
1. Introduction 

Recently, soft materials such as colloidal and granular systems have attracted 
considerable interest owing to their rich behaviours. For instance, such systems can 
be in a supercooled state under conditions of low temperatures and high densities [1]. 
Under such conditions, a characteristic time acts as a function of system parameters such 
as temperature or density, and often obeys a non-Arrhenius law, which substantially 
influences the properties of the materials. Understanding the mechanism of such 
anomalous dynamical behaviours in many-body systems is important in the field of 
statistical physics. 

Kinetically constrained spin model (KCSM) is a simple model that follows the non- 
Arrhenius law [2, 3, 4]. Thus far, it has been rigorously proved that the characteristic 
times in some kinds of KCSM on finite-dimensional lattices show super-Arrhenius type 
and Vogel-Fulcher type behaviours [5, 6] . Furthermore, in the case of a KCSM on a Bethe 
lattice, another non-Arrhenius type behaviour of a characteristic time has been found by 
Monte Carlo simulations [7]. From a static aspect, this non-Arrhenius type behaviour 
is due to a nonergodic transition corresponding to k-coie percolation; this transition 
is not a thermodynamic phase transition. Concretely, at this percolation point, the 
characteristic time diverges with a power law (non-Arrhenius law); this divergence is 
supposedly controlled by a mode-coupling equation [7]. 

The above mentioned mode-coupling equations are also believed to be related 
with the anomalous dynamical behaviours of colloidal or granular systems. A related 
conjecture is that a finite-dimensional system mimics a nonergodic transition described 
by a mode-coupling equation in a mean field sense although it is not a true nonergodic 
transition but a strong finite-size effect [1]. In the case of KCSM, it has been rigorously 
proved that the nonergodic transition observed specifically in a KCSM on a Bethe lattice 
does not occur in the model on a finite dimensional lattice although there are the strong 
finite-size effects arising from the nonergodic transition on the Bethe lattice [8, 9]. 

This leads us to consider whether the mechanisms of the appearance of such strong 
finite-size effects for different finite-dimensional systems have common features. In 
order to resolve this problem, it is necessary to find a relationship between KCSM 
and mode-coupling equations. However, the mode-coupling equation describing the 
nonergodic transition has not been derived yet for KCSM on Bethe lattices: there have 
been some related studies on the derivations of mode-coupling equations for KCSM 
[10, 11, 12, 13, 14]. 

In this paper, as a preliminary step to understand such a relationship, we attempt 
to clarify the dynamical aspect of the universality class of the nonergodic transition 
observed in Fredrickson-Andersen model (a KCSM) on a Bethe lattice. Concretely, we 
derive approximately dynamical systems from this model using a perturbation analysis 
which provides exact final states as stationary solutions above the nonergodic transition. 
We find that the universal class of the nonergodic transition cannot be captured by 
each dynamical systems even at any order by itself. Nevertheless, we find that the 
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differences between the time-dependent solutions of the derived dynamical systems and 
the results obtained by Monte Carlo (MC) simulations are systematically reduced on 
a perturbation series. Furthermore, we find that this systematic perturbation analysis 
clarifies the existence of a dynamical scaling law, which provides a implication for a 
universal relation between a size scale and a time scale near the nonergodic transition. 

2. Model 

Let us consider a regular random graph consisting oi N &N sites, each of which connects 
to c G N sites chosen randomly, where N is the set of natural numbers. Then G(c, N) 
is defined as a set of such regular random graphs. For the spin variable o"j G {—1,1} 
defined on each site i G {1, ■ ■ ■ , A^} in a graph Q G G(c, A^), the Hamiltonian we consider 
is 

1 ^ 
H{ct) = -J2^,, (1) 

1=1 
where we express cr = (<7j)iLi collectively. Here, as a preliminary step to define the 
dynamics of the system, let us consider a transition rate r(cr,FjCr) from cr to -Fjcr, 
which satisfies the detailed balance condition. Here, Fi is the spin flip operator such 
that Ficr = (cti, ■ ■ ■ , — ctj, ■ ■ ■ , ctn). Let Bi be a set of sites connected to site i. Next, we 
consider the following dynamical rule. If the number of upward spins on the sites in set 
Bi are more than or equal to /c G N, the spin on site i does not flip absolutely; otherwise, 
the spin on site i flips at a transition rate r(cr, FjCr). In other words, under this rule, the 
transition rate T(cr — )■ FjCr) from cr to FjCr is expressed by r(cr, FjCr)6(2A;— c— ^ .g^_ ctj), 
where Q{x) = 1 for a; > 0, otherwise 0. We define the situation of a spin ctj with 
Q{2k — c — XlieB "^i) = s-s 'kinetically constrained' or simply 'constrained'. The 
master equation for the probability P{cr, t) that spin configuration at time t is cr is 

N 

dtP{cr, t) = Y,[T{F.,a- -^ a)P{F,a, t) - T{a ^ F,a)P{(T, t)]. (2) 

j=i 
In this paper, we consider the case 

.H((T)-H(Fi(T)), 
r{cT,F,CT) = min(l,exp( ^ ^ ^^ ' ^^ ). (3) 

Under this constrained dynamical rule, it may be confirmed that the canonical 
distribution is a stationary distribution. Further, in equilibrium, the magnetization per 
site is mcq(T) = tanh(l/2T), and the energy density is — tanh(l/2T)/2. Therefore, there 
are no thermodynamic phase transitions in the system. In this paper, MC simulations 
are performed by the following rule. First, a site i is randomly chosen. Next, the spin 
on site i fiips with the probability T{cr — > FjCr). This step is repeated and time t = 1 is 
defined by A^ repeated steps. It is plausible that in the thermodynamic limit, this MC 
simulation is the same as the dynamics of the system described by equation (2). 

Here, we briefiy review the static aspect of a nonergodic transition in the system 
caused by the constrained dynamics, which is discussed in the previous study [7]. 
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Suppose that a spin o"j is constrained. If this constraint is permanent, we define the 
situation of a spin o"j as 'frozen'. Here, let us consider a Cayley tree, which has the 
same local structures as those of the random graph, ignoring the effects of the loop 
length O(logiV). Let g E {1, 2, ■ • ■ , (^max} be a generation of a Cayley tree where g = I 
is assigned to the root. Let us consider the probability Qg that a spin at the g-th 
generation obtained under equilibrium conditions dependent on T is frozen and upward 
without considering the state of spin at the {g — l)-th generation. From the tree structure 
of the graph, we can obtain the relation 



g,_i = F{Qg), (4) 




F{Qg) = p{T)y:{ \^ lQ^(l-g,r^-^ (5) 




where p{T) = 1/(1 + exp(— 1/T)). It should be noted that by solving recursion 
equation (5) for given values of Qg^^^, Qg for g <^ (^max becomes a solution Q{h) 
satisfying Q(T) = F{Q(T)). When c = 4 and k = 3, Q(T) is zero for sufficiently 
high temperatures. However, when the temperature is decreased, Q{T) suddenly can 
take a finite value at finite temperature T^ = 0.480898, as shown in the left-hand side 
of figure 1. This singular point is k-core percolation point, below which the system is 
nonergodic. Using the quantity Q{T), the fraction of frozen spins is described as 

Q(T)"(1 - Q(T))^-'^ 

(g'(T))"(l - Q'{T)r-\ (6) 

where Q'{T) = p{T) ^^"^^ ( ^ ~ ] Q(T)"(1 - Q{T)Y-'^-'\ In this model, it has been 

known that for 2 < /c < c, this type of nonergodic transition occurs at T = Tc with 
< Tc < oo. For k = c, Tc is zero, and for k <2,Tc is oo. The schematic phase diagram 
is shown in the right-hand side of figure 1. 

3. Simple analysis of the dynamics 

Although master equation (2) provides the complete information about the system, it 
is very difficult to extract useful information of the system from (2) because the system 
has 2^ states, which is quite a large number when N is large. To avoid this difficulty, we 
consider describing the system by iVefj G N number of variables and derive approximately 
a dynamical system closed by the set of the variables, where N^s remains finite when 
N -^ oo. For simplicity, we consider the relaxation behaviours of the system for the 
initial condition m(0) = —1 where rh{t) = ^i=i<yi/N . It should be noted that the 
spin configuration for m(0) = — 1 has no constrained spins. In the following analysis, 
we fix a graph Q G G(c, iV) for sufficiently large A^ without considering ensembles for 
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Figure 1. (left) The fraction (/) of permanently constrained spins for c = 4, fc = 3. 
(right) Schematic phase diagram. 



dtPi(o;t) 



K) 



G(c, A^). In other words, the following analysis can be applicable to almost all graphs 
Q G G{c,N) in the thermodynamic limit. 

As a first step to derive an effective dynamical system, let Pj(cr; t) be the probability 
that (jj takes a at time t and Pj((cr')|cr'; t) be the probability that the spin configuration 
on the sites in set Bi is {(y'j)jeBi provided that ai takes a' at time t. Then, we have the 
following exact evolution equation. 

0(2A:-c- J]4) 
i'eB, 

[-r,P,((T;t)P,((a;.)|a;t)+r_^P,(-or;t)P,((a;.)|-a;t)], (7) 

where we define r„ = min(l,exp(— cr/T)). Here, we assume that the value of Pi{{a'j)\(y) 
does not depend on the chosen site i. This assumption may be plausible because at least, 
in this case, inhomogeneous properties of the system arising from the effects of loops of 
the random graph Q may be negligible in the thermodynamic limit. This assumption 
corresponds to the assumption that Pj((cr^)|cr) is the same as the conditional probability 
P(((t') |(j; t) that if a site with spin variable a is randomly chosen, the spin configuration 
of its nearest neighbor sites is {cr'^). With this assumption, equation (7) is rewritten as 

dtp.it) = Y,Q{2k-c- ^4) 



K) 



j'&B„ 



[-r.p.(t)P((a;)|a;t) + r_.p_.(t)P((a;.)|-a;t)], (8) 

where B„ is a set of sites connected to a site with spin variable a and p. = J2i=i Pi{^)l^ ■ 
It is plausible that p^- is identical to ^,- 5{a — ai)/N in the thermodynamic limit N — )■ oo. 
That is, the magnetization is expressed as m{t) = J^a'^Po-i't)- 

Next, we use the approximation that the spin variables on individual sites are 
independent of each other, which is exact for equilibrium spin configurations. That is. 
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Figure 2. Time-dependent magnetization J^a'^P'^ described by (11) and to by the 
MC simulations with N = 10^ for c = 4, fc = 3. 



we rewrite P(((T^)|cr; t) as 

Further, we can obtain a simpler expression as follows. 



2k 



1:4)11 "-;(') = E 

j'€B^ jGB^ l=f 




p^,{typ+,{t) 



c-l 



(9) 



(10) 



where we define f = c — k + 1. Equations (7), (9), and (10) lead to a closed dynamical 
system in terms of p^°^ = iPa)a, 



V^=G(°nP^°^) 



:ii) 



It should be noted that stationary solutions of dynamical system (11) provide the 
exact final states above the nonergodic transition point. However, behaviours of the 
system at intermediate time scales are very different from the MC simulations, as seen 
in figure 2. This means that approximation (9) fails to capture heterogeneous spin 
configurations responsible for the dynamics in intermediate time scales. In this paper, we 
use the fourth-order Runge-Kutta method for obtaining solutions of dynamical systems 
with the time discretization dt = 10~^. 

Here, let us consider a persistent time Tp{i) of a site i, which is the time span for 
a spin to fiip at site i. Next, let us focus on the target site i in the set Cj of sites on 
which constrained spins are connected to each other. Clearly, Tp{i) strongly depends 
on the entire spin configuration of sites in set C, because the spins from a edge site 
in set Cj must be fiipped in order to fiip the target spin a,. Furthermore, the spin 
configuration on the sites in set Cj is heterogeneous in terms of the spin direction if 
it is prepared from the equilibrium spin configurations. This suggests that in order to 
detect growing relaxation times of the system with respect to the persistent time, it is 
necessary to obtain information of the heterogeneous configurations of the connected 
constrained spins. 
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4. Perturbation analysis of the dynamics 

From the resuhs in the previous section, heterogeneous spin configurations seem to 
play important roles in the growing relaxation times of the system. Therefore, in order 
to continue our analysis, we attempt to apply the information of the surrounding spin 
configurations of a target site perturbatively to effective variables by increasing the value 
of iVeff . Such a method has been previously applied to some systems and was successful 
in determining some nontrivial dynamical properties [15, 16]. 

4-1. First layer 

We attempt to apply the information of the first 'layer' of surrounding spin 
configurations of a target site to effective variables. To this end, we suppose Wi G 
{1, ■ ■ ■ , c} as the number of downward spins on sites in set -Bj. With this representation, 
site i is characterized by ((Tj,Wi). Here, let Pi{a,w;t) be the probability that {ai,Wi) 
takes (a, w) and Pjj((T, w|(t', w'; t) be the conditional probability that (ctj, Wi) takes (a, w) 
provided that {aj, Wj) takes {a', w'). Of course, a trivial relation J2a X]w=o -P«('^' ""^5 ^) — 
1 holds. Using these expressions, we have the following exact evolution equation. 

dtPi{a, w; t) = {-r^Pi{a, w; t) + r^Pi{a, w; t))e{w - f) 



c 



jeB, a' w'=f 

{Pij{cr, w + a'la', w'; t) - Pij{(T, w\(t', w'; t)), (12) 

where / is defined as c — /c + 1 and Q{x) = 1 for a; > 0, otherwise 0. 

Here, we assume that the value of Pjj(o", w|(T',t(7') does not depend on the chosen 
sites i and j if two sites i and j are chosen among the pairs of sites which have 
the same distance. This assumption may be plausible because at least, in this case, 
inhomogeneous properties of the system arising from the effects of loops of the random 
graph Q may be negligible in the thermodynamic limit. This assumption corresponds 
to the assumption that Pij{cr,w\(7',w';t) with j G Bi is the same as the conditional 
probability Pci{cT,w\a', w'; t) that after a site characterized by (o"', w') is randomly chosen 
, then one of its nearest neighbor sites, when randomly chosen, is characterized by (o", w). 

Here, let us define 

P^^[w\a ,w ]t) = p . I , — ^-^, 13 

Pci{a\a',w']t) 

where Pci{(y\(y' -.w'^t) is defined in a similar way as Pci{(y-,w\a\w'\t). In fact, we can 
obtain 

\c,(v,% (14) 

(a = -l). 
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Figure 3. Time-dependent magnetization J2a'l2'w=o'^P<^,^ described by (20) and m 
by the MC simulations with N = 10*^ for c = 4, fc = 3. 



Therefore, 

dtpa,u{t) = i-(^r^p^,^{t) + crr^Pa,«,(t))0(M - /) 

c 

+ J2Y1 ^^'P^',n,'{t)C,{w'){P:,{w + a'\a', w'; t) - P:,{w\a', w'; t)), (16) 

0-' w'=f 

where Po-.w = X]i=i -P«('^''"^)/^- ^^ is plausible that p^j^^ corresponds to J2i=i^i^ " 
ai)S{w — Wi)/N in the thermodynamic limit A^ — ?■ oo. That is, the magnetization is 
expressed as m(t) = J2a J2'w=o'^P'^,w(t)- Here, in order to obtain a closed description in 
terms of p*-^-' = {P(T,w)a,w, we use the following approximation. 

P:,{w\a',w';t) = P,\{w\a';t), (17) 

which is exact for equilibrium states. Hence, we can obtain 



P:,iw\-l;t) 

p:M + ^;t) 



TL=o^p'rAiy 

(c - w)p„^^{t) 



(19) 



Yl,=Q{c-w)paAty 

Equations (16), (17), (18) and (19) lead to a closed dynamical system in terms of p^^^ 



as follows. 



9ip«=G«(P^'^)- 



(20) 



The stationary solutions of dynamical system (20) also provide the exact final states 
above the nonergodic transition point. The dynamical behaviours of dynamical system 
(20) are closer to the MC simulations than those of dynamical system (11). However, 
as seen in figure 3, the behaviours of dynamical system (20) are also gradually deviated 
by the MC simulations if the temperature approaches the transition point T^. This 
deviation indicates that approximation (17) does not determine the behaviours at low 
temperatures. 
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4-2. Second layer 

We attempt to apply the information of the second 'layer' of the surrounding spin 
configurations of a target site to effective variables. Here, for including such effects, 
we consider the following characterization of a site. First, let /» G {l,---,c} be the 
number of downward spins on sites in set Bi. Second, let Ui G {l,---,c} be the 
number of upward-constrained spins and Vi G {1, ■ ■ ■ , c} be the number of downward 
constrained spins on the sites in set B^. Thus, a site i is characterized by {ai,li,Ui,Vi). 
With these expressions, let Pi{(7,l,u,v;t) be the probability that {(7i,li,Ui,Vi) takes 
{a,l,u,v). In addition, let Pij{a,l,u,v\a',l',u',v';t) be the conditional probability 
that {ai,li,Ui,Vi) takes {a,l,u,v) provided that {(7j,lj,Uj,Vj) takes {a' , I' , u' , v') , and 
let Pj^J2J3(((Ti, /i,-ui,t>i), (o"2, l2,U2,V2)\<73, /s, %, ^3; t) be the conditional probability that 
{ai^,li^,Ui^,Vi^) and {ai^,li^,Ui2,Vi,,) take {cri,li,ui,vi) and {(J2,l2,U2,V2), respectively, 
provided that {aig,lig,Ui^,Vi.^) takes ((T3, /3,M3,f3). Of course, a trivial relation 

Ea ELo E^=o E1=o ^i(^' ^' «. v;t) = l holds . 

As in the previous sections, we assume that the value of Pij{(r,l,u,v\a',l',u',v') 
does not depend on the chosen sites i,j if two sites i,j are chosen among the pairs 
of sites which have the same distance. In addition, we assume that the value of 
Pj^J2J3((o'i, /i,'Ui,fi), (o"2, /2, W2,f2)|o"3, /3,'U3,f3) does uot depend on the chosen sites 
ii,i2,i3 if three sites ii, ^2 and ^3 are chosen among the sets of three sites which have 
the same relationship for their distances in the order. These assumptions are plausible 
because at least, in this case, inhomogeneous properties of the system arising from 
the effects of loops of the random graph Q may be negligible in the thermodynamic 
limit. This assumption corresponds to the assumption that Pij{a,l,u,v\(7',l',u',v'), 
where j G Bi, is the same as Pci{(y,l,u,v\a',l',u',v') and to the assumption that 
Piii2i3ii^iJi,Ui,Vi),{a2,l2,U2,V2)\a3,k,U3,V3), where ^2 e 5ii,?3 e Bi,^{i3 ^ ii), is 
be the conditional probability Pc2{{o'i,h,ui,vi),{(72,l2,U2,V2)\<J3,l3,U3,V3) that a site 
characterized by (0-3, I3, U3, V3) is randomly chosen first, after which one of two connected 
sites, which is connected to the first chosen site, is randomly chosen, then among the two 
connected sites, the far site from the first chosen site is characterized by {ai,li,ui,vi), 
and the near site is characterized by {<72,h,U2,V2). 

Further, we define 

pi,a^ I / ;/ ' ,. _ Pcii(T,l,u,v\a',l',u',v') 

P^'Au,vW,l,u,v)= p^^(^^,|^.^.^.^.) • (21) 

Using this representation, we also define 

P^^'"^ {h, ui,vi\ (0-2, I2, M2, V2), {as, I3, U3, V3)) = 

-Pc2((<7l, h, Ui, Vi), ((72, /2, U2, V2)\(J3, I3, U3, V3) 



Pc2{ai, ll\{a2, I2, U2, V2), (0-3, /3, U3, V3)) 

P^'f '"' {u2, V2\(y3, h, U3, t;3)^ci(o-2, I2W3, I3, Us, V3). (22) 

In fact, we can obtain 

Pci{a,l\a',l',u,v') = -Ca^i{l',u',v'), 

c 
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C^,i{l',u',v')={ 



u' (a =!,/</) 

v' {a = -l,l<f) 

I'-v' ia = -l,l>f) 

c-l'-u' (ff = !,/>/), 

1 



(23) 



(24) 



^^02(0-1, h I ((T2, h, U2, V2), (<T3, h, M3, V3)) = -C^i,h {k, U2, V2, (T3, /g), 

C — 1 

Cai,h{k,U2,V2,(^3,h) = 
' U2 - 5(^3 - ai)0(/ - k) (ai = 1, /i < /) 

V2 - 6{as - (T,)Q{f-k) (ai = -1, k < f) 

I2-V2- 5(^3 - ai)0(/3_- /) ((Ti = -1, /l > /) 

C-I2-U2- 5(^3 - ai)0(/3 - /) (^1 = 1, /l > /). 

Using these expressions, we can write the evolution equation as follows. 

dtPa,l,u,v = {-raPa,l,u,v + r_^p_<^,;,„,„) 0(/ - /) 
c c-l' I' 

+ X] 5Z 5Z 5Z '^'y'P<y',i',u\v'C„^i{l', u', v') 

a' l'=f u'=Ov'=0 
I dLI+u' /■ \ I ii I i\ T->l,'7/i \ I rf I i\l 

[P^[ {u,v\a ,1 ,u,v)- Pj {l,u,v\a ,1 ,u,v)\ 

c c-l" I" 

+ X] X] ^^+lP+l,Z"X',«"Ca',F(^",M",^^") 
l"=fu"=Ov"=0 

a' u'=0 v'=0 

[P!^2^{u + 6{a' + l),v + 6{a'-l)\{a',f-l,u',v'),{+l,l",u",v")) 
-P^au,v\ia'J-l,u',v'),ia"J",u",v"))] 

c c-l" I" 
+ ^ ^ y^ r-lP-l,l",u",v"Ca>^l>{l", u", v") 
l"=fu"=Ov"=0 

E E E ^/"'("'' ^'1 - 1' ^"' «"' ^")ax^', «', ^', -1, n 

a' u'=0 i)'=0 

[P,'f(« + 5(a' + l),t; + 5(a'-l)|(a',/,w',^'),(-l,r«",t'")) 
-P<!f(tz,t;|(a',/,<t;'),(-l,/",«",t;"))], 

where Po-,;„t, = ^;^^^Pi(o^l^u,v\t)lN . It is plausible that Pa,iuv corresponds to 
X]i=i ^(^ ~ (^i)^{!" ~ h)^{u ~- Ui)6{v — Vi)/N in the thermodynamic limit N — )■ oo. Using 
this, the magnetization is expressed as m(t) = J2a XlLo X]m=o '^v=o ^ Pcr,iuv{'^) ■ ^"^ order 
to obtain a closed dynamical system in terms of p*^^-* = {pa,iuv)a,iuvi '^^ use the following 
approximat ions . 

^02''"' (^1. ^1 1 (^^2, ^2, M2, ^^2), (0-3, ^3, M3, ^^3)) 

= P'r\h,u,,v,\a2,l2'^ f), (26) 

(27) 



(25) 



^/,(T / 



,/,(T / 



,i (m, t;|(7 , / , M , t; ) = P,i (m, t;|a , / ^ /), 
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Figure 4. Time-dependent magnetization J2a X]/=o ^1=0 '^v=o '^Pcy.iuv described by 
(32) and m by the MC simulations with N = 10'^ for c = 4, fc = 3. 



which are exact for equihbrium states. In fact, we can obtain the concrete expression 
of PcP^'"(m, v\a', I'^f) as follows. 



P^i< ''^(m, v\ + lj' > f) = {c-l- u)p„^iuvl ^ (c - / - u)p^^iuv, (2 



Kf,uv 



H^^^'^iu, V\ -l,l'>f) = {l- V)p,,luj 5^ (/ - V)p,, 
P^^^ '''{U,V\ + !,/'</)= Vpa^luvl ^ Vp^^luv, 



luvi 



(29) 
(30) 
(31) 



where T.i>f,uv = EL/ ELo E^=o and Ei</,™ = E/=o E^=o E^=o- That is, equations 
(25) — (31) lead to a closed dynamical system in terms of p'-^^ as follows. 

V')=G(2)(P^'^)- (32) 

The stationary solutions of dynamical system (32) also provide exact final states of 
the system above the nonergodic transition point. Further, the solutions of dynamical 
system (32) are closer to the MC simulations than those of previous two descriptions. 
However, similar to previous descriptions, the solutions are gradually deviated from the 
MC simulations for lower temperatures, as seen in figure 4. 

5. On a perturbation series 

5.1. Systematic improvement of the solutions 

We consider the manner of changes in the solutions of the derived dynamical systems on 
a perturbation series. It is noteworthy that the differences between the MC simulations 
and the time-dependent solutions of the derived dynamical systems are systematically 
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Figure 5. The magnetization described by each dynamical systems and rh by the MC 
simulations at T = 1. TV = 10^, c = 4, fc = 3. 



reduced with the increase in the number n of p^"'\ as observed in the left-hand side of 
figure 5. 

The MC simulations show that under the initial condition that all spins are 



downward, the magnetization behaves as rh{t) 



m. 



eq 



(T) ~ exp(— t/r) in a long time 



limit, as shown in the right-hand side of figure 5. As shown in figure 6, a rough estimation 
of r indicates the behavior r ~ e~^, where C — 3 and e = {T/T^ — 1). In fact, the value 
around C — 3 has been already confirmed for the persistent time in a previous study [7]. 
Let us consider the relaxation time T„(r) of the dynamical systems defined as 
|_\(")|-i^ where A*-"-* is the maximum eigenvalue, except for trivial zero, of matrix 
A^(") obtained by linearizing G'"^ at a stationary solution p^ = limi_5.oo p''"^^) with 
m(0) = —1. Tn{T) does not show the power-law behaviour e~^, as seen in figure 6 
although when e is slightly below 1, Tn(T) behaves as e~^" with (n < Cn+i < C? and 
Tn{Tc) is finite. Thus, higher order analyses will be needed for capturing the power- law 
behavior e""". 



5.2. Systematic construction of higher order perturbations 

In this section, we consider a perturbation analysis that is of a higher order than those 
discussed thus far. Let us reiterate the way to characterize a site discussed in section 
4.1. In the analysis, in order to define the effective state of a target site i, we use 
information of sites in set Bi. We regard the way of such a characterization of site i 
as the first order characterization. Next, let us reiterate the way to characterize a site 
discussed in section 4.2. In the analysis, in order to define the effective state of a target 
site i, we use information of sites in set Bi and set Bj with j G B^. In other words, 
site i is characterized by the number of sites in set -Bj, where each site in set Bi is 
characterized by the information of the first order characterization. We regard the way 
of such a characterization of site i as the second order characterization. 

Using this second order characterization, we can characterize the sites in set B^ 
except for the information of the branch directed to site i. Next, site i is characterized 
by the number of sites in set 5j, where each site in set B^ is characterized by the above 
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Figure 6. e-dependence of each relaxation time. MC simulations are performed with 
N = 10^, c = 4,fc = 3. 



characterization. This procedure defines third order characterization of site i. In the 
same way, we can define arbitrary n-th order characterization iteratively. 

Further, essentially the same approximation as those used in the first order and 
second order perturbation analysis can be applicable to n-th order perturbation analysis, 
which provides exact final states of the system as stationary solutions above the 
transition point. With this procedure, in principle, we can compute the relaxation time 
Tn{T) described by the dynamical system for all n orders. In this perturbation series, the 
formal 'small' parameter can be regarded as the distance between which sites are used 
for defining the effective states of a target site. It should be noted that the concrete 
value of n does not have any significance other than in the sequence, and oo-order 
perturbation analysis may provide the original Master equation in the thermodynamic 
limit by the definition. 



5. 3. A dynamical scaling law in maximum eigenvalues at T^ 

The systematic improvement of the solutions described by the derived dynamical 
systems on a perturbation series and the systematic construction of a higher order 
perturbation series motivate us to consider some scaling relation between order n in 
the perturbation series and relaxation time (the inverse of the maximum eigenvalue ) 
Tn(T). The first question raised here is how large n causes the divergence of r„(Tc) 
as lim„_j,„^ r„(Tc) — )■ oo. On the basis of the consideration that the size of connected 
constrained sites with the heterogeneous spin configuration resulting in the maximum 
relaxation time can be infinite for equilibrium spin configurations at T = T^, we 
can expect Uc = oo. Explicitly, using a function F dependent on (c, fc), we express 
Tn(^c) = F{n, c, k) where F{n, c, fc) — t- cxo with n — )■ cxd. 

As shown in the figure 6, since tq seems to be out of the scaling region even if 
some scaling relation exists, we focus on ti(Tc) and T2(Tc) for various parameters (c, k). 
Let us remind that this perturbation analysis includes the effects of larger connected 
constrained spins if order n is increased. Therefore, it may be plausible that accessible 
correlation sizes by the perturbation analysis with order n are increased as order n is 
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Figure 7. ^^(rc) as a function of (7?,(c, fc; Q)n") ^ where we set 7?.(4, 3; a) = 1, 
7?,(5,4;q) = 1 and a = 1. Fitting parameters are TZ{5,3; 1) — 0.51, 7^(6,3; 1) = 0.37, 
7^(7, 3; 1) = 0.3, 7e(6, 4; 1) = 0.48 and n{7, 4; 1) = 0.35. (•, •) in the figure is (c, fc). 



increased. On the basis of this consideration, first, we assume that TZ{c, k; a)"!^" expresses 
accessible correlation sizes by the perturbation analysis with order n, where 7^ and a are 
certain constants. It should be noted that the determination of a needs the information 
about the nature of correlation sizes in the system, which is discussed in section 5.4. 
Furthermore, we also assume that r„ has power-law forms in the accessible correlation 
sizes in the system. Therefore, if we set a value of a, the above assumptions lead to the 
exponents of the power-law forms. That is, we assume the following form: 

TniTc) = J■a(7^(c, t, «)n", c, k), (33) 

where 7^(c, k; a) are fitting parameters, and J^a{x, c, k) = Cx^* . Surprisingly, as seen in 
the figure 7, we can find that there are values of 7^ for arbitrary values of a such that 
J^a{x, c, k) is independent of c within the numerical analysis. This result indicates that 
the assumption for the power-law dependences of r„ on order n is plausible. Concretely, 
az*(fc;a) ~ z'{k) {z'{3) c^ 2.15), z'{4) ~ 2.4). In the following, we present some 
conjectures related to the nature of correlation sizes in the system. 



5.4- Conjectures arising from dynamical scaling law (33) 

The law (33) implies that there exists a characteristic size TV obeying the dynamical 
scaling law t{T) ~ N'{TY where t{T) has been already defined through the relaxation 
of the magnetization. In fact, we have already known a candidate of A/'(T), which is 
called the minimal rearrangement size defined as the minimal number of flipped spins 
on the surrounding sites of a target spin in order to flip the target spin [17, 18]. The 
previous study has captured that the size scale shows the power-law behaviour such as 
^-'^ic,k) ng^j. |;]2g nonergodic transition in FA model [18]. 

Here, we mention to the meanings of parameter a. Rough numerical simulations 
indicate that we need set a = amrs ~ 1.1 in order to identify z* as z. In this context. 
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c^mrs can be regarded as the quantity connecting the perturbation order n to accessible 
minimal rearrangement size by the perturbation analysis with order n. In other words, 
if we obtain a exact value of ttmrs such that amrs-2 = z' in a. certain case of (c, k) = (c', k'), 
the values of z{= ^'/amrs) in any cases (c, k) 7^ (c', k') can be derived by the perturbation 
analysis presented above, because the analysis provides z' in any cases of {c,k), in 
principle. 

The result for z* leads to the conjecture that z also does not depend on the value of 
c with the same value of k. This means that the universality classes of nonergodic 
transitions in the kinetically constrained spin model are classified by the value of 
constraint parameter k, and the quantity characterizing the universality is not ( or 
u but z, where ( and u depend on the value of c. Actually, we have performed the MC 
simulations in order to confirm the above conjecture. Although we have found signs 
supporting the above conjecture, we could not obtain plausible results due to the finite 
size effects and the limitation of the maximum step of time. It is an important future 
study to confirm the conjecture by MC simulations. 

6. Concluding remarks 

In this study, we have constructed a systematic perturbation analysis for the dynamics of 
FA model on a Bethe lattice. This systematic perturbation analysis clarifies the existence 
of a dynamical scaling law, which provides a implication for a universal relation between 
a size scale and a time scale near the nonergodic transition. 

Here, we discuss the relevance of our results to the previous studies. Actually it has 
been conjectured that the persistent time found by MC simulations can be described by 
a mode-coupling equation [7]. This statement is not inconsistent with Uc = 00 because 
mode-coupling equations are 00-dimensional differential equations. In addition, a fact 
supporting the relevance of the model to a mode-coupling equation has been made in 
the literature of the analysis of the minimum size rearrangement [17]. In addition to the 
previous results, the results obtained in this paper provide another plausible conjecture 
for the properties of the nonergodic transition. That is, the nonergodic transitions have 
a weak universality. This means that critical exponents (, v for time and size depend 
on the value of c, but the dynamical critical exponent z does not depend on the value 
of c with the same value of k. A similar statement has been mentioned in the previous 
study for the dynamical transition in j9-spin glass model on the Bethe lattice [19]. The 
study states that z does not depend on the quantity 7 related to the connectivity of 
the graph, which plays the similar role to that of c in this paper. However, the study 
does not mention to the dependence of z on the value of j9, which may play the similar 
role to that of k in this paper. The results in this paper indicate that z depends on the 
value of p. The confirmation of this conjecture for the dynamical transition in p-spin 
glass model is an important future study. 

Another aspect of this weak universality appears in the comparison to the case of a 
ferromagnetic Ising model on a Bethe lattice. In the previous study, the critical exponent 
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can be obtained using the dynamical system closed by finite number of effective variables, 
which is derived by the similar approximation method to that of this paper [15]. 
Therefore, the universality class of nonergodic transition in the KCSM is quantitatively 
and qualitatively different from that of ferromagnet-paramagnet transition in some spin 
models including, at least, Ising model. Here, the word 'qualitatively' means that the 
differences are located in not only the value of critical exponents and the existence of the 
dynamical system closed by finite number of variables capturing the critical exponent. 
Finally, we discuss the relationship between nonergodic transitions discussed in this 
paper and the related phase transitions in other systems [20]. In fact, a decimation 
dynamics of a random graph in the thermodynamic limit exhibits a saddle-node 
bifurcation at the k-core percolation point [21, 22]. Furthermore, a random-field 
Ising model with zero-temperature Glauber dynamics in the thermodynamic limit 
at a spinodal transition has been reported to correspond to the k-core percolation, 
which is also a saddle-node bifurcation on Bethe lattices [24, 25]. It should be noted 
that dynamical behaviours near these transitions corresponding to the saddle-node 
bifurcation are extremely different from the dynamics near the nonergodic transition 
in FA model although those are k-core percolations from the static viewpoint. This 
difference may be related to the existence of a no-passing property in the system [23] . 
That is, FA model does not have a no-passing property whereas the other systems 
do. This consideration leads to the conjecture that the nonergodic transition in other 
kinetically constrained spin models such as Kob- Andersen models belong to universality 
classes that are different from that of the saddle-node bifurcation; Such models do 
not have the no-passing property [3, 9]. Finally, we mention a previous study that 
suggests a relationship between the static properties of jamming transition and the k- 
core percolation in a previous study [26, 27]. Since such a system showing the jamming 
transition do not have the no-passing property, finding a relationship between the 
macroscopic dynamical behaviors near the jamming transition and KCSM is also an 
interesting topic for future studies. 
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